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TITLE OF THE INVENTION 

OPTIMUM DESIGN METHOD AND APPARATUS, AND 
PROGRAM FOR THE SAME 

BACKGROUND OF THE INVENTION 

Field of tfre in v ention 

[0001] The present invention relates to a solution of an 
optimization problem of design parameters, and more 
particularly to an automatic design technique for 
optimization of shapes including topology of structural 
members . 

D escri p tion of the Related Art 

[0002] The term "structural topology optimum design" 
means a problem for deciding the topology and shape 
dimensions of a structural member, which are optimum under 
given conditions. Hereinafter, the topology and shape 
dimensions of a structural member will be referred to as 
"design argument functions", and the above decision problem 
will be referred to as a "design-variable-function 
optimization problem". The reason why the term "argument 
functions" is used resides in that the topology and shape 
dimensions are given as functions of a three-dimensional 
space . 
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[0003] In the design-variable-function optimization 
problem, an optimization problem of a state argument 
function must be solved for a value of each design argument 
function. From this point of view, the structural topology 
5 optimum design can be regarded as a double structure 

optimization problem that contains an optimization problem 
of a state argument function on the inner side and an 
optimization problem of a design argument function on the 
outer side. As to the optimization problem of a state 
10 argument function on the inner side, the concept of dividing 

a space into a finite number of elements is employed based 
on accumulation of known techniques . 

[0004] For a problem taking strain energy of a structural 
member as an evaluation generic function, in particular, a 
15 finite element method is generally employed as a solution 

technique for that problem. A direct method with respect to 
a linear equation is employed for solution of the finite 
element method. 

[0005] On the other hand, regarding the optimization 
20 problem of a design argument function, the following three 

kinds of methods have been primarily provided (see Reference 
1, S. Bulman, J. Sienz, E. Hinton: "Comparisons between 
algorithms for structural topology optimization using a 
series of benchmark studies". Computers and Structures, 79, 
25 pp. 1203-1218 (2001), or Reference 2, Y-L. Hsu, M-S. Hsu, C-T. 
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Chen: "Interpreting results from topology optimization using 
density contours". Computers and Structures, 79, pp. 1049- 
1058 (2001)): 

1. Evolutionary method (referred to as "E method" below) 
5 2. Homogenization method (referred to as "H method" below) 

3. Material distribution method (referred to as "MD 
method " below ) 

[0006] According to the E method, each of partial spaces 
resulting from dividing a space is called a cell, and 

10 creation and disappearance of cells are repeated in 

accordance with proper rules . A structural member is given 
as a set of cells that exist finally. A definite structural 
member is obtained by allowing only two states represented 
by whether the cell exists or not. Also, because the E 

15 method does not employ differential information of the 

evaluation generic function and can avoid trapping to a 
local optimum solution, the E method is effective in the 
case of the evaluation generic function having multiple 
peaks . 

20 [0007] Reference 3 (Japanese Patent Laid-Open No. 2001- 

134628) proposes an optimization design device for a 
skeleton structural member, which employs one kind of the E 
method, i.e., a genetic search method. The proposed 
optimization design device comprises an approximate 

25 optimization computing device using approximation formulae 
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for discrete design variable data such as section dimensions 
of a skeleton structural member, and a close optimization 
computing device using the design variable data- Those two 
computing devices are combined with each other to be adapted 
5 for an actual design problem containing a large number of 

design variables. 

[0008] The H method enables sensitivity analysis to be 
employed by assuming a finer structure for a structural 
element positioned in each of divided partial spaces and 

10 introducing a new design argument function having a 

continuous value. Here, the term "sensitivity analysis" 
means an analytical technique utilizing differential 
information of the evaluation generic function with respect 
to the design argument function. If the sensitivity 

15 analysis is enabled, it becomes possible to employ an 

iterative method, such as a gradient method. As a result, a 
computation time required for search of at least a local 
optimum solution can be greatly cut in comparison with a 
round-robin technique, such as the E method. (See Reference 

20 4, Hiroshi Yamakawa, "Saiteklka Dezain (Optimization 

Design)", Computation Dynamics and CAE Series 9, Baifukan 
Co. , Ltd. (1996) . ) 

[0009] The MD method is a method for expressing changes 
in topology and shape dimensions of a structural member by 
25 assigning, to each element, a real number in the range of 0 
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to 1 representing the presence probability of the structural 
member. The MD method is similar to the H method in that 
the sensitivity analysis is enabled by replacing discrete 
information, which indicates whether the structural member 
5 is present or not, with a continuous value representing the 

presence probability of the structural member. However, 
because the number of parameters is smaller in the MD method 
than in the H method, the MD method is advantageous in that 
modeling is easier and the application range is wider. 

10 [0010] Reference 5, Fujii, Suzuki and Otsubo, "Structure 

Phase Optimization Using Voxel Finite Element Method" , 
Transactions of JSCES, Paper No. 200000010(2000), describes 
a phase optimization technique for a structure based on the 
MD method. This technique has features given below. 

15 [0011] (1) Because of using a voxel finite element method 

( in which a space is divided at equal intervals ) , element 
rigidity matrices for all elements are the same. 
Accordingly, by computing the element rigidity matrix once 
in advance, it is available in subsequent computations. 

20 Furthermore, since elements are regularly arranged, there is 

no need of storing node number information for each element. 
[0012] (2) A pre-processed conjugate gradient method and 
an element-by-element method are combined with each other to 
solve large-scaled simultaneous linear equations. As a 

25 result, a solution can be obtained without setting up an 
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entire rigidity matrix, and therefore the memory capacity 
required for processing can be reduced. 

[0013] (3) A homogenization method requires 6 design 
variables (in the three-dimensional case) for one element. 
5 Further, each time the design variable is changed, the 

element rigidity matrix must be computed again. On the 
other hand, by employing a density method expressing the 
presence probability of the structural member as a density 
ratio, only one design variable is required for one element. 
10 In addition, changes of the design variable do not affect 

the element rigidity matrix. 

[0014] The above-mentioned known techniques, however, 
have problems given below. 

[0015] Generally, a structure optimization problem is 
15 formulated as a double optimization problem containing an 

optimization problem of a state variable vector per 
iterative step for an optimization problem of a design 
variable vector. Assuming the optimization problem of a 
state variable vector to be optimization on the outer side 
20 and the optimization problem of a design variable vector to 

be optimization on the inner side, the optimization on the 
inner side is a problem of determining the state variable 
vector with the design variable vector being a parameter, 
i.e., the design variable vector being fixed. This problem 
25 is usually called structural analysis and can be solved 
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based on the finite element method by using the solution 
technique for a linear equation. 

[0016] However, if a structure changes and a structural 
member does not exist in a certain area any more, for 
5 example, if a structural member is holed, the design 

variable of the corresponding element becomes 0 and the 
Young's modulus of that element also eventually becomes 0. 
This leads to a result that the problem cannot be solved by 
the direct method because a coefficient matrix of the linear 
10 equation is not fully filled with elements and an inverse 

matrix does not exist. 

[0017] For that reason, most of the known techniques 
employ a method of using a value of not exactly 0, i.e., 
replacing 0 with a small value close to 0, when the design 

15 variable vector becomes 0. This method is similarly applied 

to the voxel finite element method employed in Reference 5. 
[0018] Setting the design variable of the element to a 
small value, however, corresponds to the case in which there 
exists a thin film or a weak member from a physical point of 

20 view. Stated another way, an area where no materials are 

present is not exactly expressed in the known techniques. 
[0019] Further, as seen from the figure showing the 
result in Reference 5, a structural element area not 
contributing to total strain energy, i.e., a floating island 

25 area or a projection area, is left in th result. Such a 
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structural element area should be finally removed, but a 
method for removing that area is not explained in Reference 
5. 

[0020] In Fig. 16, by way of example, when elements 512, 
5 513 and 514 are neither subjected to any weight nor 

supported, those elements have no contributions from the 
standpoint of strength. In other words, those elements are 
ones not contributing to total strain energy and hence 
should be removed in accordance with a design demand for 
10 minimum total weight of the structural members. However, it 

is difficult to find those elements so long as the solution 
technique is based only on simple geometrical information, 
i.e., on connection of characteristic functions. 

15 SUMMARY OF THE INVENTION 

[0021] It is an object of the present invention to 
provide an optimum design method capable of executing 
structure optimum design, including topology, without 
20 utilizing the experiential rules that depend on individual 

problems . 

[0022] Another object of the present invention is to 
provide an optimum design method capable of eliminating a 
component corresponding to a structural element that does to 
25 contribute to an evaluation function, when a variable for 
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optimizing a value of the evaluation function is determined. 
[0023] According to one aspect, the present invention 
which achieves those objects relates to an optimum design 
method comprising a first solution determining unit of 
5 solving an optimization problem of a first evaluation 

function for a state variable vector with a design variable 
vector being as a parameter; and a second solution 
determining unit of solving an optimization problem of a 
second evaluation function for the design variable vector 

10 and the state variable vector obtained in the first solution 

determining unit , wherein the second solution determining 
unit includes a gradient vector computing unit of computing 
a gradient vector of the second evaluation function for the 
design variable vector; a first coefficient computing unit 

15 of computing a first coefficient based on a value of a norm 

of the gradient vector; a search vector computing unit of 
computing a search vector based on the first coefficient; a 
second coefficient computing unit of computing a second 
coefficient; and a design variable vector updating unit of 

20 updating the design variable vector based on the second 

coefficient, the second coefficient computing unit including 
the first solution determining unit, the first solution 
determining unit being executed as an iterative method based 
on the gradient vector, and the state variable vector being 

25 not initialized during iteration. 
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[0024] According to one aspect, the present invention 
which achieves those objects relates to an optimum design 
apparatus comprising first solution determining means for 
solving an optimization problem of a first evaluation 
5 function for a state variable vector with a design variable 

vector being as a parameter; and second solution determining 
means for solving an optimization problem of a second 
evaluation function for the design variable vector and the 
state variable vector obtained in the first solution 

10 determining means, wherein the second solution determining 

means includes gradient vector computing means for computing 
a gradient vector of the second evaluation function for the 
design variable vector; first coefficient computing means 
for computing a first coefficient based on a value of a norm 

15 of the gradient vector; search vector computing means for 

computing a search vector based on the first coefficient; 
second coefficient computing means for computing a second 
coefficient; and design variable vector updating means for 
updating the design variable vector based on the second 

20 coefficient, the second coefficient computing means 

including the first solution determining means, the first 
solution determining means being executed as an iterative 
method based on the gradient vector, and the state variable 
vector being not initialized during iteration. 

25 [0025] According to still another aspect, the present 
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invention which achieves those objects relates to a 
computer-readable optimum design program, the program 
comprising codes for causing a computer to perform an 
optimum design method comprising a first solution 
5 determining step of solving an optimization problem of a 

first evaluation function for a state variable vector with a 
design variable vector being as a parameter; and a second 
solution determining step of solving an optimization problem 
of a second evaluation function for the design variable 

10 vector and the state variable vector obtained in the first 

solution determining step, wherein the second solution 
determining step includes a gradient vector computing step 
of computing a gradient vector of the second evaluation 
function for the design variable vector; a first coefficient 

15 computing step of computing a first coefficient based on a 

value of a norm of the gradient vector; a search vector 
computing step of computing a search vector based on the 
first coefficient; a second coefficient computing step of 
computing a second coefficient; and a design variable vector 

20 updating step of updating the design variable vector based 

on the second coefficient, the second coefficient computing 
step including the first solution determining step, the 
first solution determining step being executed as an 
iterative method based on the gradient vector, and the state 

25 variable vector being not initialized during iteration. 



[0026] According to yet another aspect, the present 
invention which achieves those objects relates to an optimum 
design method comprising a first solution determining step 
of solving an optimization problem of a first evaluation 
function for a state variable vector with a design variable 
vector being as a parameter; and a second solution 
determining step of solving an optimization problem of a 
second evaluation function for the design variable vector 
and the state variable vector obtained in the first solution 
determining step, wherein the second solution determining 
step includes a design variable vector updating step of 
updating the design variable vector in sequence, the design 
variable vector updating step including a minimum point 
searching step of making search from a start point to obtain 
a minimum point; and a terminal point evaluating step of 
deciding an optimum point based on a value of the second 
evaluation function at the minimum point and a value of the 
second evaluation function at an end point. 

[0027] According to another aspect, the present invention 
which achieves those objects relates to an optimum design 
apparatus comprising first solution determining means for 
solving an optimization problem of a first evaluation 
function for a state variable vector with a design variable 
vector being as a parameter; and second solution determining 
means for solving an optimization problem of a second 
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evaluation function for the design variable vector and the 
state variable vector obtained in the first solution 
determining means , wherein the second solution determining 
means includes design variable vector updating means for 
5 updating the design variable vector in sequence, the design 

variable vector updating means including minimum point 
searching means for making search from a start point to 
obtain a minimum point; and terminal point evaluating means 
for deciding an optimum point based on a value of the second 

10 evaluation function at the minimum point and a value of the 

second evaluation function at an end point . 
[0028] According to a further aspect, the present 
invention which achieves those objects relates to a 
computer- readable optimum design program, the program 

15 comprising codes for causing a computer to perform an 

optimum design method comprising a first solution 
determining step of solving an optimization problem of a 
first evaluation function for a state variable vector with a 
design variable vector being as a parameter; and a second 

20 solution determining step of solving an optimization problem 

of a second evaluation function for the design variable 
vector and the state variable vector obtained in the first 
solution determining step, wherein the second solution 
determining step includes a design variable vector updating 

25 step of updating the design variable vector in sequence, the 



design variable vector updating step including a minimum 
point searching step of making search from a start point to 
obtain a minimum point; and a terminal point evaluating step 
of deciding an optimum point based on a value of the second 
evaluation function at the minimum point and a value of the 
second evaluation function at an end point. 
[0029] According to a further aspect, the present 
invention which achieves those objects relates to an optimum 
design method comprising a first solution determining step 
of solving an optimization problem of a first evaluation 
function for a state variable vector with a design variable 
vector being as a parameter; a second solution determining 
step of solving an optimization problem of a second 
evaluation function for the design variable vector and the 
state variable vector obtained in the first solution 
determining step; and an erasing step of erasing, from the 
design variable vector, a component corresponding to a 
structural element that does not contributes to the second 
evaluation function . 

[0030] According to a further aspect, the present 
invention which achieves those objects relates to an optimum 
design apparatus comprising first solution determining means 
for solving an optimization problem of a first evaluation 
function for a state variable vector with a design variable 
vector being as a parameter; second solution determining 
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means for solving an optimization problem of a second 
evaluation function for the design variable vector and the 
state variable vector obtained in the first solution 
determining means; and erasing means for erasing, from the 
5 design variable vector, a component corresponding to a 

structural element that does not contributes to the second 
evaluation function . 

[0031] According to a further aspect, the present 
invention which achieves these objects relates to a 

10 computer-readable optimum design program, the program 

comprising codes for causing a computer to perform an 
optimum design method comprising a first solution 
determining step of solving an optimization problem of a 
first evaluation function for a state variable vector with a 

15 design variable vector being as a parameter; a second 

solution determining step of solving an optimization problem 
of a second evaluation function for the design variable 
vector and the state variable vector obtained in the first 
solution determining step; and an erasing step of erasing, 

20 from the design variable vector, a component corresponding 

to a structural element that does not contributes to the 
second evaluation function. 

[0032] Other objectives and advantages besides those 
discussed above shall be apparent to those skilled in the 
25 art from the description of a preferred embodiment of the 
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invention which follows. In the description, reference is 
made to accompanying drawings, which form a part thereof, 
and which illustrate an example of the invention. Such 
example, however, is not exhaustive of the various 
5 embodiments of the invention, and therefore reference is 

made to the claims which follow the description for 
determining the scope of the invention. 

BRIEF DESCRIPTION OF THE DRAWINGS 

10 

[0033] Fig. 1 is a flowchart showing a processing 
sequence of a second solution determining step. 
[0034] Fig. 2 is a block diagram of an apparatus 
according to one embodiment . 
15 [0035] Fig. 3 is a flowchart showing a processing 

sequence of a line searching step. 

[0036] Fig. 4 is a flowchart showing a processing 

sequence of a search range narrowing step. 

[0037] Fig. 5 is a flowchart showing a processing 
20 sequence of a minimum point searching step. 

[0038] Fig. 6 is a flowchart showing a processing 

sequence of a first vector correcting step. 

[0039] Fig. 7 is a flowchart showing a processing 

sequence of a second vector correcting step. 
25 [0040] Fig. 8 is a flowchart showing a processing 
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sequence of a coefficient B computing step. 
[0041] Fig. 9 is a flowchart showing a processing 
sequence of a normal vector computing step. 
[0042] Fig. 10 is a flowchart showing a processing 
5 sequence of a first solution determining step. 

[0043] Fig. 11 is an illustration for explaining 
practical setting of a problem. 

[0044] Fig. 12 is an illustration showing a calculation 
result in the embodiment . 
10 [0045] Fig. 13 is a flowchart showing a processing 

sequence for removing an undesired element. 

[0046] Fig. 14 is a flowchart showing another processing 
sequence for removing an undesired element. 

[0047] Fig. 15 is an illustration showing a calculation 
15 result after removing the undesired element. 

[0048] Fig. 16 is an illustration for explaining 
undesired elements. 



DE SCRIPTION QF TflE PREFERRED jEMBQPIMSNTS 

20 

[0049] A preferable one embodiment according to the 
present invention will be described in detail below with 
reference to the accompanying drawings. 
( First Embodiment ) 
25 [0050] Fig. 2 shows the configuration of an apparatus 
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according to one embodiment- As shown in Fig. 2, the 
apparatus comprises a CPU 201, a display 202, an input unit 
203, a primary storage 204, a secondary storage 205, a 
communication unit 206, and a bus line 207. 
5 [0051] Various processing functions executed in this 

embodiment is stored as programs in the secondary storage 
205 beforehand. Any of the stored programs is loaded into 
the primary storage 204 and then executed in response to a 
command input via the input unit 203 or the communication 
10 unit 206. 

[0052] For the sake of explanation given below, a target 
problem is formulated as follows. 

[0053] Assuming that an argument function is expressed by 
a finite dimension vector with formulation based on the 

15 finite element method, an evaluation generic function for an 

argument function is given as an evaluation function for a 
variable vector. The following description is made on an 
assumption that the argument function is expressed by a 
finite dimension vector. 

20 [0054] A state variable vector x and a design variable 

vector f are written as column vectors, respectively, as 
given below; 

(Eq. 1) x = (x 1# x 2 xJ T 

(Eq. 2) f = (f 1# f 2 fJ T 

25 wherein T represents transposition, x represents an m- 
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dimension vector, and f represents an n-dimension vector. 
[0055] Boundary conditions of x and f are assumed to be B x 
and B 2 , respectively, as follows: 

(Eq. 3) B^x) = 0 

5 (Eq. 4) B 2 (f) = 0 

[0056] Evaluation functions for the state variable vector 
x and the design variable vector f are assumed to be L x and 
L 2 , respectively, as follows; 

(Eq. 5) L x = Mx.-f) 

10 (Eq. 6) L 2 = L 2 (f,x) 

wherein h 1 is a function with x being a variable vector and f 
being a parameter, and L 2 is a function with x and f being 
both variable vectors . 

[0057] In a general optimization problem, a number K x of 
15 equality restraint conditions and a number K 2 of inequality 

restraint conditions are applied in many cases. Therefore, 
the j-th equality restraint condition and the k-th 
inequality restraint condition for the design variable are 
assumed to be Q., and R k , respectively, as follows: 
20 (Eq. 7) Q 3 (f,x) = 0 

(Eq. 8) R*(f,x) * 0 

[0058] From the definition described above, the optimum 
design problem is given as a problem of determining a 
solution, which satisfies the following equation, under the 
25 restraint conditions (Eq. 4), (Eq. 5), (Eq. 6), (Eq. 7) and 
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(Eq. 8): 

(Eq. 9) min[L 2 ] 
[0059] Additionally, the state variable x is obtained as 
a solution satisfying the following equation under the 
5 restraint condition (Eq. 3): 

(Eq. 10) min[L x ] 
[0060] Processing executed in this embodiment will be 
described below on the basis of the above formulation. 
[0061] For simplicity of the explanation, it is assumed 
10 in the following description that the equality restraint 

conditions of (Eq. 7) is given by one condition as follows; 

(Eq. 11) Q(f ) = 2Cjf t ( j) + c 0 = 0 
wherein c., (j = 0, 1, — , n) is a constant. 

[0062] Also, the inequality restraint conditions of (Eq. 
15 8) are assumed to be given as the following ranges of value; 

(Eq. 12) Rj (f) = f t (j) - f 3min * 0 

V<f ) = f jma x - f t (j) ^ 0 
wherein j = 1, 2,..., n. 

[0063] Accordingly, the number of the inequality 
20 restraint conditions is twice the number of the design 

variables . 

[0064] Flag vectors a* and a" having the same dimensions 
as those of the design variable vector are defined as 
follows: 
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a* = 0, if 0 s: R/(f) 
a + = 1 , otherwise 

and 

a" = 0, if 0 ^ R^f) 
5 a" = 1 , otherwise 

[0065] Also, it is assumed that a set of j, at which a + 
and a" are each 1 , is Ai and a set of other indices is A 2 . 
[0066] Fig* 1 is a flowchart showing a processing 
sequence of a second solution determining step. Referring 

10 to Fig. 1, in step S101, the CPU reads various factors of a 

system to be simulated. The various factors may be read as 
data entered from the input unit 203 or the communication 
unit 206, or may be given by reading data that has been 
stored as a file in the secondary storage 205 beforehand. 

15 The various factors of the system include initial values of 

x and f , the boundary conditions B lr B 2 , the evaluation 
functions L 1# L 2 , and the restraint conditions Q d , R k . Based 
on such information of the various factors, the program 
secures a required variable area in the primary storage 204 

20 and sets values. 

[0067] In step S102, the number of times t of iterations 
is initialized to 1. The number of times t of iterations 
means the number of times at which a first iteration process 
has been executed in the second solution determining step. 

25 The second solution determining step is brought to an end if 
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the number of times of iterations of the processing from 
step S103 to step Sill reaches a setting value T. 
[0068] A gradient vector computing step is executed in 
step S103. In the gradient vector computing step, the CPU 
5 computes a gradient vector g t of a second evaluation function 

for the design variable vector. The gradient vector g t is 
computed in accordance with the following equation: 

(Eq. 13) g t = (dL a /3f lf dl* 2 /d 2 dL 2 /df n ) T 

[0069] Whether 3L 2 /3f j can be analytically solved or not 

10 depends on cases. If an analytical solution is present, (Eq. 

10) is obtained by utilizing a computing process for the 
analytical solution. If otherwise, the gradient vector can 
be computed by utilizing the technique called automatic 
differentiation. The automatic differentiation technique is 

15 known from, for example. Reference 6 (Koichi Kubota and 

Masao Iri, "Arugorizumu no Jldobibun to Oyou ("Automatic 
Differentiation and Applications of Algorithm" , Modern 
Nonlinear Science Series, Corona Publishing Co., Ltd., 
(1998) ) . 

20 [0070] In step S104, a first vector correcting step is 

executed for g t obtained from (Eq. 13). The first vector 
correcting step will be described later with reference to 
Fig. 6. 

[0071] A convergence determining step is executed in step 
25 S105. In this step, the square of a norm of g t is computed. 
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and if the computed value becomes 0 or smaller than a very 
small value e close to 0, the processing is brought to an 
end. The square of the norm is computed in accordance with 
the following equation; 
5 (Eq. 14) |g t || 2 = g t T g t 

[0072] A first coefficient computing step is executed in 
step S106. In the first coefficient computing step, a 
coefficient p is computed in accordance with (Eq. 15) given 
below: 

10 (Eq. 15) P = ||g t || 2 / ||g t J| 2 

wherein p = 1 holds at t = 1 . 

[0073] A search vector computing step is executed in step 
S107. In this step, a search vector p t at the t-th iteration 
is computed in accordance with (Eq. 16) given below: 
15 (Eq. 16) p t +- -g t _ x + Pp t . 2 

[0074] In step S108, a first vector correcting step is 
executed for p t obtained from (Eq. 16). This step will be 
described later with reference to Fig. 6. 

[0075] A line searching step is executed in step S109. 
20 In this step, a coefficient 6 for updating the design 

variable is decided. This step will be described later with 
reference to Fig. 3. 

[0076] A design variable vector updating step is executed 
in step S110. In this step, the design variable vector is 
25 updated as follows: 
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(Eq. 17) f t+1 = f t + 6p t 

[0077] In step Sill, a second vector correcting step is 
executed for the design variable vector. This step will be 
described later with reference to Fig. 7. 

[0078] In step S112, t is updated to t+1, and if t after 
the update exceeds the preset number of times of iterations, 
the processing is brought to an end. 

[0079] The foregoing is the flow of the second solution 
determining step. Individual processing steps in this step 
will be described in more detail below. 

[0080] The line searching step will be described below 
with reference to Fig. 3. In this step, the current design 
variable value f t and the search vector p t at the current 
design variable value are given as arguments. It is also 
assumed that the evaluation function and the restraint 
condition are given as appropriate. 

[0081] In step S301, a maximum step size A is computed 
using the following equation; 

(Eq. 18) A = MIN [ ( f f(j) - f^j) ) / p t (j) ] 

wherein ff(j) = f )max . if p t (j) > 0 
= f jmln . if p t (j) < 0 
[0082] Further, various variables are initialized as 
follows ; 

(Eq. 19) gl = (3.0 - 5.0 1/2 ) / 2 
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g2 = (5.0 1/2 - 1) / 2 

dA = A/T m 

al = 0.0 

a2 = al + dA 

5 wherein T m is a positive integer used for dividing the search 

range into appropriate zones and is set to a value in the 
range of, e.g., 10 to 100. 

[0083] In step S302, a search range narrowing step is 
executed. This step will be described later with reference 

10 to Fig. 4. 

[0084] In step S303, 6 set in the search range narrowing 
step is compared with the maximum step size A. If they are 
equal to each other, the processing is brought to an end. 
If otherwise, the CPU proceeds to step S304. 

15 [0085] In step S304, a minimum point searching step is 

executed. This step will be described later with reference 
to Fig. 5. The line searching step is thereby brought an 
end. 

[0086] The search range narrowing step will be described 
20 below with reference to Fig. 4. 

[0087] In step S401, values of the second evaluation 
functions are computed as follows; 

(Eq. 20) L21 = L2(x(f t + alp t ), f t + alp t ) 

(Eq. 21) L22 = L2(x(f t + a2p t ), f t + a2p t ) 

25 wherein x(f t + alp t ) and x(f t + a2p t ) are obtained by 
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executing a first solution determining step that will be 
described later with reference to Fig. 10. 

[0088] In step S402, L21 and L22 are compared with each 
other. If L21 is larger than L22, the processing is brought 
5 to an end. If otherwise, the CPU proceeds to step S403. 

[0089] In step S403, t is initialized to 0 and \i{t) is 
initialized to 0 . 

[0090] In step S404, \i(t) is updated as follows: 
(Eq. 22) \i(t) = |ui(t-l) + dA 

10 [0091] In step S405, L21 and L22 are updated as follows; 

(Eq. 23) L21 = L22 

(Eq. 24) L22 = L 2 (x(f x + fx(t)p t ), f x + jx(t)p t ) 

wherein x(f x + fi(t)p t ) is obtained by executing the first 
solution determining step that will be described later with 

15 reference to Fig. 12. 

[0092] In step S406, L21 and L22 are compared with each 
other. If L21 is not larger than L22, the CPU proceeds to 
step S408, and if otherwise, it proceeds to step S407. 
[0093] In step S407, t is updated to t+1, and thereafter 

20 t is compared with T m . If t is larger than T m , the CPU 

proceeds to step S409, and if otherwise, it proceeds to step 
S404. 

[0094] In step S408, computations of the following 
equations are executed: 
25 (Eq. 25) al = n(t) - 2dA 
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a2 = \x(t) 

[0095] In step S409, jJ-(t) and A are compared with each 
other. If \x{t) is larger than A, A is substituted for 6. 
If otherwise, an arbitrary value other than A is substituted 
5 for 5, whereby the search range narrowing step is brought to 

an end. 

[0096] The minimum point searching step will be described 
below with reference to Fig. 5. 

[0097] In step S501, a computation of the following 
10 equation is executed: 

(Eq. 26) da = a2 - al 

[0098] If da is not larger than a preset small positive 
real number, the CPU proceeds to step S504, and if otherwise, 
it proceeds to step S502. 
15 [0099] In step S502, computations of the following 

equations are executed: 

(Eq. 27) vl = al + gl da 

v2 = al + g2 da 
[0100] Further, the search range is narrowed as follows: 
20 (Eq. 28) a2 = v2 , if L21 < L22 

al = vl , otherwise 
[0101] This narrowing is a range reduction based on the 
golden section method. 

[0102] In step S503, values of the second evaluation 
25 function are computed based on the following equations, and 
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the CPU proceeds to step S501; 

(Eq. 29) L21 = L2(x(f x + al p t ), f x + al p t ) 

L22 = L2(x(f x + a2 p t ) , f x + a2 p t ) 
wherein x(f x + alp t ) and x(f x + a2p t ) are obtained by 
5 executing the first solution determining step that will be 

described later with reference to Fig. 10. 

[0103] In step S504, a middle point ac of the search zone 
and a value L2c of the second evaluation function at the 
middle point are computed as follows : 
10 (Eq. 30) ac = (al + a2)/2 

(Eq. 31) L2c = L2(x(f x + ac p t ) , f x + ac p t ) 

[0104] In step S505, a value of 6 is decided as follows: 
(Eq. 32) 5 = ac, if L2c = MIN[L21 ,L22 ,L2c,L2e] 

= al, else if L21 = MIN[L21 ,L22 ,L2c,L2e] 
15 = a2, else if L22 = MIN [L21,L22,L2c, L2e ] 

= A, else if L2e = MIN[L21 ,L22 ,L2c,L2e] 
= 0 , otherwise 
[0105] In (Eq. 32), MIN[ ] is a function of returning a 
minimum value among the arguments . The minimum point 
20 searching step is thereby brought to an end. 

[0106] The first vector correcting step executed in each 
of step S104 and step S108 will be described below with 
reference to Fig. 6. 

[0107] A second coefficient computing step is executed in 
25 step S601. In this step, values of a unit normal vector NN, 
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a length DST of a vertical line from the origin to a 
hyperplane, and a coefficient B are computed- The second 
coefficient computing step will be described in detail later 
with reference to Fig. 8. 
5 [0108] A first vector projecting step is executed in step 

S602. In this step, a vector X is corrected in accordance 
with an equation given below by using the unit normal vector 
NN, the length DST of the vertical line from the origin to 
the hyperplane, and the coefficient B, which have been 

10 computed in step S601: 

(Eq. 33) X(j) = X(j) - B NN(j), if j 6 A, 

= 0 , otherwise 
[0109] The first vector correcting step is thereby 
brought to an end. 

15 [0110] The second vector correcting step executed in step 

Sill will be described with reference to Fig. 7. 
[0111] The second coefficient computing step is executed 
in step S701. This step will be described in detail later 
with reference to Fig. 8. 

20 [0112] A second vector projecting step is executed in 

step S702. In this step, the vector X is corrected in 
accordance with an equation given below by using the unit 
normal vector NN, the length DST of the vertical line from 
the origin to the hyperplane, and the coefficient B, which 

25 have been computed in step S701: 
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(Eq. 34) X(j) = X(j) - B NN(j), if j 6 A, 

= f,max. if a* = 0 

= iW if a" = 0 

[0113] The second vector correcting step is thereby 
5 brought to an end. 

[0114] The second coefficient computing step will be 
described below with reference to Fig. 8. 

[0115] In step S801, a normal vector computing step is 
executed to obtain NN and DST. This step will be described 
10 in detail later with reference to Fig. 9. 

[0116] In step S802 # the coefficient B is computed using 
the following equation: 

(Eq. 35) B= ^ u( j ) NN ( j ) - DST 

[0117] The second coefficient computing step is thereby 

1 5 brought to an end . 

[0118] The normal vector computing step executed in step 
S801 will be described below with reference to Fig. 9. In 
this step, the unit normal vector NN and the length DST of 
the vertical line from the origin to the hyperplane are 

20 computed. 

[0119] In step S901, the variables expressed by the 
following equations are computed; 

(Eq. 36) DST = -c(0) ' 

wherein c(0)' is calculated in accordance with an equation 
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given below; 

(Eq. 37) c(0)'=c(0)+ ^c(k)f t (k) 

k€=A 2 

(Eq. 38) D = ( Jc(k) 2 r 1/2 

(Eq. 39) J] c( J )f t ( 3> + c( 0 ) ' = 0 

5 [0120] DST in Eq. 36 is corrected in accordance with an 

equation given below; 

(Eq. 40) DST = DST / D 

[0121] DST given by above (Eq. 40) provides the length of 
the vertical line from the origin to the hyperplane. 
10 [0122] In step S902, a normal vector computing step is 

executed. 

[0123] First, a normal vector N not normalized is 
computed in accordance with (Eq. 41) and (Eq. 42) given 
below: 

15 (Eq. 41) N( j) = ]1 c ( k ) / c ( 3)- if jG A, 

(Eq. 42) N(j) = f t (j), otherwise 
[0124] By using (Eq. 41), a unit normal vector NN for a 
partial space is computed in accordance with an equation 
given below; 

20 (Eq. 43)NN( j) = N( j)(2 11 c ( k f / c ( i ) 2 )" 1/2 , if jE A, 
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[0125] For example, when the total sum of elements of the 
design variable vector is 1 as the equality restraint 
condition, components of the unit normal vector is 
simplified as follows because of c(0) = -1, c(j) = 1, j = 1, 
5 2 t ... , n i 

(Eq. 44) NN(J) = | A, | " 1/2 , if jE A, 

wherein | A x | is the number of elements of A x . 
[0126] The normal vector computing step is thereby 
brought to an end. 

10 [0127] The first solution determining step will be 

described below with reference to Fig. 10. 
[0128] If a structure changes, there is a possibility 
that an element having a characteristic function value of 0 
appears. In the case of carrying out a structural analysis 

15 for such a structure, a coefficient matrix becomes not 

regular and the problem cannot be solved by the direct 
method. For that reason, one of the following three methods 
is employed in the first solution determining step: 

(1) The characteristic function value being 0 is replaced 
20 with a small value close to 0, and then simultaneous linear 

equations are solved by the direct method. 

(2) Simultaneous linear equations are reconstructed using 
only the elements having values above 0, and then the direct 
method is employed to solve the problem. 

25 (3) The problem is solved by an iterative method using a 
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gradient vector, which does not require a computation of an 
inverse matrix. This iterative method includes a steepest 
descent method, a conjugate gradient method, etc. 
[0129] The first solution determining step using the 
5 conjugate gradient method will be described below. 

[0130] In step S1001, the CPU reads various factors of a 
system to be simulated. The various factors include an 
initial value of x, a value of the design parameter f , the 
boundary condition B 1# and the evaluation function L x . Based 
10 on such information of the various factors, the program 

secures a required variable area in the primary storage 204 
and sets values. 

[0131] In step S1002, the number of times t of iterations 
is initialized to 1. The number of times t of iterations 

15 means the number of times at which a first iteration process 

has been executed in the first solution determining step. 
The first solution determining step is brought to an end if 
the number of times of iterations of the processing from 
step S1003 to step S1010 reaches a setting value. 

20 [0132] In step S1003, the CPU computes a gradient vector 

g 1<t of a first evaluation function for the state variable 
vector. The gradient vector g 1#t is computed in accordance 
with the following equation: 

(Eq. 45) g lt = (d'L 1 /dx 1 . dl^/dx^..., dh 1 /dx n ) T 

25 [0133] Whether dL 1 /3x j can be analytically solved or not 
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depends on cases. If an analytical solution is present, (Eq. 
45) is obtained by utilizing a computing process for the 
analytical solution. If otherwise, the gradient vector can 
be computed by utilizing the technique called automatic 
5 differentiation. The automatic differentiation technique is 

known from, for example, the above-cited Reference 6, etc. 
[0134] In step S1004, the square of a norm of g lmt is 
computed, and if the computed value becomes 0 or a value 
very close to 0, the processing is brought to an end. The 

10 square of the norm is computed in accordance with the 

following equation ; 

(Eq. 46) ||g 1>t | 2 = g./g^ 

[0135] In step S1005, a coefficient P is computed in 
accordance with (Eq. 47) given below; 

15 (Eq. 47) 0 = ||g lft || 2 / (g^f 

wherein (3 = 1 holds at t = 1 . 

[0136] In step S1006, a search direction vector p t is 
computed based on (Eq. 48) given below; 
(Eq. 48) p 1#t = -g liM + pPi.t-i 

20 [0137] In step S1007, a line searching step is executed. 

This step is to decide a coefficient a for updating the 
design variable. When the structural analysis problem is 
solved by using the finite element model, a value of a is 
given by the following equation: 

25 (Eq. 49) a = P x . t T gi. t / (Pi.jAPi.t) 
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[0138] In step S108, the state variable vector is updated 
as follows : 

(Eq. 50) x t+1 = x t + a p l t 

[0139] In step S1009, t is updated to t+1, and if t after 
5 the update. exceeds the preset number of times of iterations, 

the processing is brought to an end. 

[0140] The first solution determining step is thereby 
brought to an end. 
( Practical Example ) 

10 [0141] In this practical example, the present invention 

is applied to optimum shape automatic design for a 
cantilevered beam that is subjected to a weight at an 
arbitrary position. For the sake of simplicity, the 
following description is made on an assumption being limited 

15 to a problem of plane strain. 

[0142] As shown in Fig. 11, a design area in which a 
structural member is able to exist is a rectangle 1102. 
According to the finite element method, the design area is 
divided at equal intervals into numbers n y , n x of partial 

20 areas, respectively, in the vertical (length) and horizontal 

(width) directions. The divided partial areas are called 
cells and numbered such that a lower left cell and an upper 
right cell are expressed respectively by (1,1), (n y ,n x ). 
Lattices points are called nodes and likewise numbered such 

25 that a lower left node and an upper right node are expressed 
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respectively by (1,1) , (n y +l,n x +l). 

[0143] In Fig. 11, numeral 1101 denotes a support member, 
and 1103 denotes a weight vector. Each cell (j,k) has a 
corresponding characteristic function value f(j,k). Here, 
5 the term "characteristic function value" means a variable 

taking a positive real number value in the range of 0 to 1 , 
which represents the existence probability of the structural 
member in the cell (j,k), and it is an element of the design 
variable vector f , given below, in the present invention: 

10 (Eq. 51) f = (f(l,l), f(l,2),..., f(n y ,n x )) T 

[0144] Similarly, each node (j,k) has a corresponding 
horizontal displacement u(j,k) and vertical displacement 
v(j,k). These displacements are real numbers taking 
arbitrary values and are elements of the state variable 

15 vector U, given below, in the present invention: 

(Eq. 52) U = (u(l,l), v(l,l), u(l,2), v(l,2) 

u(n y +l,n x +l) , v(n y +l,n x +l) ) T 
[0145] Assuming a rigidity matrix to be A and a weight 
vector to be b, as a well-known result based on the finite 

20 element method, the state variable vector U is given as a 

solution of the optimum problem of an evaluation generic 
function expressed by the following equation: 

(Eq. 53) L 2 m (1/2) U T AU - b T U 
[0146] More specifically, it is known that the state 

25 variable vector U is given as a solution of a linear 
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equation given below: 

(Eg. 54) AU - b = 0 
[0147] Further, g x t computed in the first solution 
determining step is given as the left side of (Eq. 54). 
Because (Eq. 54) is a set of simultaneous linear equations, 
it is usually solved by the direct method based on an 
inverse matrix. In the structure optimum problem, however, 
the rigidity matrix A is a function of the design variable 
vector f . Therefore, if f contains an element taking a 
value of 0, the rigidity matrix A is not fully filled with 
elements and becomes not regular. This results in that the 
simultaneous linear equations cannot be solved by the direct 
method. For that reason, the present invention employs the 
conjugate gradient method. 

[0148] An evaluation function L 2 for the design variable 
vector f is defined by total strain energy as follows: 

(Eq. 55) L 2 = (1/2)U T AU 
[0149] Then, g t required in the second solution 
determining step can be computed in accordance with the 
following equation ; 

(Eq. 56) g t = dL 2 /3f(e) = -(1/2) U e T A e U e 
wherein U e , A e are respectively a displacement vector at a 
node corresponding to an element e and an element rigidity 
matrix corresponding to the displacement vector. 
[0150] Fig. 12 shows a computation result obtained in 
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this example. In Fig. 12, a black area represents an area 
in which the structural member exists. 

[0151] In this example, the aspect ratio of the design 
area for the structural member is 2 : 1 (length/width) and 
5 its analytical solution is, as known, a combined structure 

of beams arranged at ± 45 degrees relative to the horizontal 
direction. 

[0152] It is understood that the computation result shown 
in Fig. 12 is well matched with the known analytical 

10 solution. 

[0153] Thus, according to the embodiment described above, 
the structure optimum design including topology can be 
executed without utilizing the experiential rules that 
depend on individual problems . 

15 [0154] An undesired element removing process will be 

described below. Fig. 13 is a flowchart showing a 
processing sequence for removing an undesired element. 
Referring to Fig. 13, in step S1301, the CPU reads various 
factors of a system to be simulated. The various factors 

20 may be read as data entered from the input unit 203 or the 

communication unit 206, or may be given by reading data that 
has been stored as a file in the secondary storage 205 
beforehand. The various factors of the system include 
initial values of x and f , the boundary conditions B lf B 2 , 

25 the evaluation functions L 1# L 2 , and the restraint conditions 
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Q j# R k . Based on such information of the various factors, 
the program secures a required variable area in the primary 
storage 204 and sets values, 

[0155] In step S1302, a solution of the optimization 
5 problem under the inequality restraint conditions formulated 

as (Eq. 4) to (Eq. 9) is determined. While the above- 
described method can be used in this solution determining 
step, the undesired element removing process described below 
is also applicable to the case of employing any of the 
10 following methods: 

( 1 ) sequential linear plan method 

(2) feasible direction method 

(3) gradient projection method 

(4) general contraction gradient method 
15 (5) optimum reference method 

[0156] Subsequent steps S1303 to S1308 constitute the 
undesired element removing process. 

[0157] First, in step S1303, an element s is initialized 
to 1. 

20 [0158] In step S1304, it is checked whether f(s) is 0. 

If so, the CPU proceeds to step S1308, and if otherwise, it 
proceeds to step S1305. 

[0159] In step S1305, the value of f(s) is changed as 
follows : 

25 (Eq. 57) f(s) <— f(s) + 6 
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wherein 6 is an arbitrary real number other than 0 , at which 

f ( s ) takes a value not smaller than 0 , but not larger than 1 . 

[0160] In step S1306, sensitivity dL 2 /df(s) of a second 

evaluation function for f(s) is computed as follows: 
5 3L 2 /df = d/df [ (1/2)U T AU] 

[0161] From a rigidity equation of AU = b, the following 

equation is obtained: 

(Eq. 58) dh 2 /df = d/df [ (l/2)b T U] = (l/2)b T dU/df 

[0162] On the other hand, assuming the weight b to be 
10 constant without depending on f , partial differentiation of 

both sides of the rigidity equation of AU = b by f leads to: 
(Eq. 59) au/df = -A-'dA/dfU 

because of (dA/df )U + A(dU/df) = 0. 

[0163] By putting (Eq. 59) in (Eq. 58), 
15 (Eq. 60) dL 2 /df = -(l/2)U T (dA/df )U 

is obtained from the relation of U = A _1 b. 

[0164] The above equation can be rewritten for a 

characteristic function of the element s as follows; 
(Eq. 61) dL 2 /df(s) = -(l/2)U s T (dA s /df (s))U s 
20 wherein U s is a vector constituted by a displacement on the 

node belonging to the element s . and A s is an element 

rigidity matrix corresponding to U s . 

[0165] Whether dA e /df (e) can be analytically solved or not 
depends on cases. If an analytical solution is present, a 
25 computing process for the analytical solution can be 
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obtained. In the case of A s = f(s)A s ', for example, dL 2 /df(s) 
is given by the following equation. 

(Eq. 62) dL 2 /df(s) = - ( 1/2 )U S T A S ' U s 
[0166] If dA/af(s) cannot be analytically solved, it is 
5 possible to execute the computation by utilizing the above- 

mentioned technique called automatic differentiation. 
[0167] Here, the sensitivity expressed by (Eq. 62) is 
given as 0 or a small value close to 0 for all the elements 
as a result of the solution determining step executed in 
10 step S1302. 

[0168] A description will be made of whether the value of 
(Eq. 62) is changed or not when the characteristic function 
value f(s) of the element s is changed in accordance with 
(Eq. 57). 

15 [0169] First, when the element s contributes to the value 

of strain energy, i.e., when (1/2)U 8 T A S U S is not 0, A s changes 
with a change of f(s), and therefore the strain energy 
(l/2)U s T A g U s of the element s also changes. Thus, U s obtained 
as the solution of the structural analysis problem takes 

20 different values U s ' between before and after the change of 

f(s). Accordingly, (Eq. 62) takes a negative value as 
follows : 

(Eq. 63) dL 2 /3f(s) = - ( 1/2 ) U s ' T A S » U s ■ * - ( 1/2 )U S T A S ' U s = 0 
[0170] On the other hand, when the element s does not 
25 contribute to the value of strain energy, i.e., when 



- 42 - 



(1/2)U S T A S U S is 0, U s becomes 0 because of A s being a positive 
definite matrix, and therefore (1/2)U S T A S U S does not change 
even with a change of f(s). Accordingly, U s obtained as the 
solution of the structural analysis problem takes the same 
value between before and after the change of f(s), whereby 
(Eq. 62) becomes 0. 

[0171] Stated another way, whether the relevant 
structural element is necessary or not can be determined by 
checking whether an absolute value of the sensitivity 
resulting upon a change in the value of f(s) is 0 or a small 
value close to 0 . 

[0172] In step S1307, whether the sensitivity dh 2 /df(s) is 
0 is checked. If so, f(s) is updated to 0, and if otherwise, 
f(s) is not updated. 

[0173] In step S1308, s is updated to s+1. Then, if 
updated s exceeds n, the processing is brought an end, and 
if otherwise, the CPU proceeds to step S1304. 
[0174] In the method described above, the undesired 
element removing step is executed as a post-process. 
However, when the second solution determining step is 
executed as an iterative process, the undesired element 
removing step may be executed in the iteration loop. The 
processing in such a case will be described below with 
reference to Fig. 14. The conjugate gradient method is 
employed, by way for example, as the iterative process in 
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the following description. 

[0175] In step S1401, the CPU executes initialization. 
More specifically, the design variable vector f is set to a 
predetermined value, which written as f 0 . 
5 [0176] In step S1402, t is set to 1. 

[0177] In step S1403, a value which the gradient vector g t 
of the second evaluation function for the design variable 
vector f takes at f = f 0 is computed as follows: 
(Eq. 64) g t ■ 3L 2 /df 

10 = (dL 2 /3f(l), dL 2 /df(2) dL 2 /af(n)) T 

[0178] In step S1404, g t is corrected so as to satisfy the 
equality restraint condition and the inequality restraint 
condition . 

[0179] In step S1405, it is checked whether a norm 
15 computed in accordance with the following equation exceeds a 

preset value. If the norm exceeds the preset value, the 

processing is brought to an end, and if otherwise, the CPU 

proceeds to step S1406; 

(Eq. 65) ||g t || = (g t T g t r 1/2 

20 [0180] In step S1406. p* defined by an equation given 

below is computed: 

(Eq. 66) p = ||g t || / 

wherein {3 = 0 holds at t = 1 . 

[0181] In step S1407, the search vector p t is computed in 
25 accordance with the following equation: 
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(Eq. 67) p t = pp t "g t 

[0182] In step S1408, p t is corrected so as to satisfy the 
equality restraint condition and the inequality restraint 
condition. 

5 [0183] In step S1409, line search is executed along p t to 

find f that minimizes the second evaluation function. The 
found f is set as f t . 

[0184] In step S1410, f t is corrected so as to satisfy the 
equality restraint condition and the inequality restraint 
10 condition. 

[0185] In step S1411, the undesired element removing step 
is executed. 

[0186] In step S1412, t is updated to t + 1, and if t after 
the update exceeds a preset value, the processing is brought 

15 to an end. If otherwise, the CPU proceeds to step S1403. 

[0187] Additionally, the processing executed in steps 
S1404, S1408 and S1410 can be made by using a technique 
described as the gradient projection method in the above- 
cited Reference 4. 

20 (Practical Example) 

[0188] In this practical example, the above -described 
embodiment is applied to optimum shape automatic design for 
a cantilevered beam that is subjected to a weight at an 
arbitrary position. For the sake of simplicity, the 

25 following description is made on an assumption being limited 
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to a problem of plane strain. 

[0189] The evaluation function L 2 for the design variable 
vector f is defined by total strain energy in accordance 
with (Eq. 55) mentioned above: 
5 (Eq. 55) L 2 = (1/2)U T AU 

[0190] Then, f minimizing (Eq. 55) is determined. This 
kind of problem is usually subjected to, as the equality 
restraint condition, a restriction on the total weight being 
constant, namely; 

j-riy k = n x 

10 (Eq. 68) ^ ( j,k)= constant 

j-l k = l 

and to, as the inequality restraint condition, a restriction 
on the value range which the characteristic function value 
can take, namely: 

(Eq. 69) 0 s f(j,k) <; 1 
15 [0191] As described above, that optimization problem 

under the inequality restraint condition can be solved by 
the known solution determining method. 

[0192] The undesired element removing step in this 
practical example will be described below. 
20 [0193] In particular, the computation of sensitivity 

executed in step S1306 of Fig. 13 can be made by using the 
following equation; 

(Eq. 70) 3L 2 /af(j,k) = - ( 1/2 )U 3 . k T A d , k U jJC 
wherein ll, ik is an element displacement vector having, as a 
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component, a displacement on the node belonging to the 
element (j,k), and A., k is an element rigidity matrix 
corresponding to the element displacement vector. 
[0194] Subsequently, a floating island element can be 
5 removed in accordance with the method described above with 

reference to Fig. 13. 

[0195] Fig. 15 shows the shape of the structural member 
after the undesired element removing step has been executed 
on the shape of the structural member shown in Fig. 12. It 

10 can be confirmed that floating island areas and projections 

appearing in Fig. 12 are eliminated in Fig. 15. In this 
example, the undesired element is determined and removed 
only when the absolute value of the sensitivity is exactly 
equal to 0 . A computing time required for the undesired 

15 element removing step is about 40 seconds when the CPU is 

Pentium III (933 MHz). 

[0196] Incidentally, the present invention is applicable 
to not only an apparatus constituted by a single unit, but 
also to a system constituted by a plurality of units. 

20 Further, the present invention may be implemented by 

supplying, to an apparatus or a system, a storage medium 
storing program codes of software for realizing the 
functions of the above -described embodiment, and by causing 
a computer in the apparatus or the system to read and 

25 execute the program codes stored in the storage medium so 
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that the apparatus or the system achieves the desired 
functions . 

[0197] In addition to the case of directly realizing the 
functions of the above- described embodiment by the computer 
5 in the apparatus or the system executing the program codes 

stored in the storage medium, the present invention also 
involves the case of realizing the functions of the above- 
described embodiment by causing an OS (Operating System) , 
etc. operating on the computer to execute the processing in 

10 accordance with instructions of the program codes. 

[0198] In such a case, a storage medium storing the 
program codes constitutes the present invention. 
[0199] Although the present invention has been described 
in its preferred form with a certain degree of particularity, 

15 many apparently widely different embodiments of the 

invention can be made without departing from the spirit and 
the scope thereof. It is to be understood that the 
invention is not limited to the specific embodiments thereof 
except as defined in the appended claims. 



